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Field-theoretical method is proposed, that yields analytic expression for lateral 
pressure distribution across hydrophobic core of a bilayer lipid membrane. Lipid 
molecule, repelling entropicaly surrounding neighbors, is modeled as a flexible 
string embedded between confining walls in the self-consistent harmonic poten- 
tial. At room temperature the pressure profile is expressed solely via lowest 
energy eigenfunctions of the operator of the string energy density. Analytical 
results compare favorably with molecular dynamics simulations. The theory is 
suited to study protein-lipid interaction. 
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Development of microscopic theory of fluid lipid membrane is a nontrivial problem. One 
of the fundamental difficulties lies in the elongated flexible structure of the constituting lipid 
molecules that results in their large conformational entropy. Neighboring lipids essentially re- 
duces the manifold of accessible conformations of the molecule's hydrocarbon chains in a bilayer 
membrane as compared with a free lipid molecule. This reduction, in turn, leads to effective 
entropic repulsion between lipid molecules of non-local nature. Entropic interaction^ between 
hydrocarbon chains causes unusual for common fluids inhomogeneous lateral pressure distribu- 
tion in thermodynamic equilibrium state. The lateral pressure profile strongly influences the 
functioning of proteins embedded in the membrane^, ^, for instance by shifting their activation 
barriers^^,^^. Due to a nanometer scale of membrane thickness the pressure proflle is difficult 
to measure. Several attempts to attack the present problem theoretically have already been 
made. Some of them involved molecular dynamics simulations of a lipid bilayer*,"^. The others 
exploited numerical methods (e.g. Monte Carlo) with mean-field approach, introducing lat- 
eral pressure profile as a substitute for inter- molecular entropic repulsion^, In the present 
work we use path-integral technique in combination with adiabatic invariants theorem^ and 
find lateral pressure profile in analytical form. Due to its generality our method permits, in 
principle, analytical derivation of the lateral pressure profile for a great variety of membrane 
constituents and lipid variations (head group type, chain length, etc.) by appropriate change 
of the imposed boundary conditions and membrane spontaneous curvature. Hence, our theory, 
besides fundamental interest, may be used as a tool for theoretical prediction of the properties 
of a wide variety of lipid membranes. A sketch of lipid bilayer membrane is given in Fig.l. 
Hydrocarbon chains form hydrophobic part in lipid molecules. Collisions between chains to- 
gether with excluded volume effect give rise to entropic repulsion between lipid molecules. This 
repulsion is responsible for inhomogeneous lateral pressure in the hydrophobic core of the mem- 
brane. The ensemble of chains is described within mean- field approximation using microscopic 
model of a single hydrocarbon chain in effective potential, see 1. The chain is modeled as a 
flexible string (see Fig. 2) that possesses bending rigidity, Kf. The effective potential, which 
allows for entropic repulsion between neighboring chains in each monolayer is expressed via 
lateral pressure distribution H(z) at the depth z multiplied by area itB?{z) swept by the chain 
in {x,y} -plane:t/e// — 7rHi?^(z) (compare e.g.^). Here the chain deviation from the z-axis 
is = + Ry. The energy functional of a single hydrocarbon chain Et is chosen in the 
following form: 

Here L is chain length; z is coordinate along the chain axis; Ris vector in {x, y} plane character- 
izing chain deviation from z-axis. In the functional (1) the bending energy term contains second 
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derivative over z-coordinate rather than over the contour length of the chain. The approxima- 
tion made is vaUd provided that chain deviations from the z-axis are small with respect to chain 
length L (applicability of this condition at room temperature is justified below). We rewrite 
the energy functional Ei using self-adjoint operator, if , that is obtained after integrating by 
parts the expression in (1): 

^* = E f Ri{z)HRi{z)dz -,11= lif/^ +7rn(^). (2) 

The boundary conditions are: i) chain angle is fixed in the head group region; ii) zero torque 
is applied at the free chain end; iii) no total force is applied to the chain at head group and at 
the free chain end. These conditions are expressed as: -R'(O) = R"'{0) = R"{L) = R"'{L) = 0. 
An arbitrary chain conformation described by " shape- functions" Rx,y{z) can be expanded over 
orthonormal set of the eigenfunctions, Rn{z), of operator H: 

Rx = ^CnRn{z). (3) 

Since x- and y-deviations are independent we shall retain for simplicity only one index, e.g. a;. 
Then, the energy functional Et of hydrocarbon chains given in (2) takes equivalent form: 



Et= f Y.^-^nRl{^)dz, 
Jo „ 
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where En are the eigenvalues of the operator H in (2). 

To find the distribution of lateral pressure along 2-axis we shall use the following formal 
trick^. Imagine that we press the membrane laterally at depth z and measure its response. This 
would give us local lateral pressure. Mathematically this procedure can be performed using the 
adiabatic invariants theorem^. Consider the energy (refenhe) as a sum of potential energies 
of harmonic oscillators labeled by index n and coordinate z that have rigidities kn = 'i.En 
and vibration amplitudes = CnRnizY"^ . Then, the lateral pressure is a sum of individual 
contributions produced by the oscillators defined above. We calculate each contribution as a 
thermodynamically averaged force exerted by a one-dimensional oscillator on the surrounding 
"walls". Assume that the oscillator is centered at x = and the walls are positioned at points 
±(i(t), where d{t) — d — vt and v 0. Here d is z-indepcndent radius of a cylinder, which 
on average hosts a single chain in a flat bilayer, see Fig. 2. The rate of the energy change 
dcnz/dt of n-th oscillator at depth z is related to the force, fn,z, exerted by the oscillator on 
the wal\s:d€nz / dt = 2vfn,z- Alternatively, the rate of the energy change can be calculated using 
the method of adiabatic invariants^: 

di/dd 

denz dt = V—— (5) 

dl/den,z 



where time-independent adiabatic invariant is :/ — Jpdx. Here p ~ ■\/2m(e„_z — knx'^/2) is 
momentum, x is coordinate of the oscillator (mass ni of the oscillator does not enter the final 
expression). Integration in I gives area inside the closed orbit in phase-space {p,x} crossing 
at points ±d(t) the x-axis and at ±p(d{t)) the p-axis. Equating the above two expressions for 
dcnz/dt and taking limit w ^ we find the force: 

/{\/knenz\/^ ^ d'^kn/2enz Ml. ic) ^ i 
v2arcsin(d-\/fc„/2e„j;) (."J 
0, rf2fc„/2e„^ > 1. 

The force produced by n-th oscillator at depth z is found as ci result of averaging over Cn'- 

< fn{z) >= r U^nz{Cn))P{Cn)dCn, (7) 
J —oo 

where e„2(C„) = £^„C^i?^(z) and the distribution function equals: 
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PiCn) ^ (E^/nkBT)'/' oxp 1^--^ E„C^Rl{z)dzj. (8) 
Summation over n permits to find the total force at depth z: 



(9) 



Obtained formula is remarkable: though probabilities for "oscillators" n with different 2;'s 
are distributed independently according to formula (8), nevertheless, factors Rn{z) in formula 
(9) keep track of the fact that "oscillators" arc different degrees of freedom of one and the same 
vibrating string. Actually, only the oscillators that reach the walls at ±d contribute to the force 
f{z). This bounds energy e from bclow:e > (Pkn/2. Hence, the lower limit of integration in 
(9) increases fast with index n, because bending part of total energy in (1) results in power- law 
dependence: fc„ = 2_E„ cx n^. Simultaneously, there is an upper limit of integration in (9) 
~ R'^{z)kBT imposed by the Boltzmann factor. Hence, non-zero integration interval exists 
only for the first few terms in (9). Here we limit ourselves to the terms n = 0, 1 and find the 
thermodynamically averaged force: 

</W>-E^^nW- (10) 
n 

As the membrane is laterally isotropic, this force can be considered as a radial force in the 
{x — y} plane. Then, the pressure distribution is obtained dividing the radial force in (10) by 
half the perimeter nd : 

nW«E^^n(^)- (11) 

n 

This is the central result of our paper, since, we had shown that the lateral pressure profile 
that enters the chain energy functional Et in Eq. (1) is expressed merely via the eigenfunctions 
of the energy density-operator H. Using this result we find the corresponding self-consistent 
solution for Il{z) analytically at room temperature. From the estimates for the chain bending 
rigidity: Kf ~ kbTlp, where Ip is chain persistence length, it proves that the separation between 
the eigenvalues £"„ of the operator H is large with respect to k^T. This justifies approximate 
projection of the unknown eigenfunctions Rn{z) on the reduced basis formed by the lowest 
energy eigenfunctions of the free-chain operator Hq = K f /2{d'^ /dz): 

Ro = aiYo + biYi- Ri = azFo + 62^1, (12) 

where: 

Yo{z) = l/vT, Yn>i{z) = ^/2/L{cos{Xnz) + cos{XnL) cosh.{X„z) / cosh(A„L)}, 
and the eigenvalues are given by: 

Wo = 0, Wn>i = Kf/{2L^){TT/4)\4n + 3)\ 

where A„>iL = 37r/4 -|- nn. The coefficients a^, bi satisfy the conditions of orthonormality of 
the functions Rn{z)- 

af + bf = 1, aibi + 0262 = 0. 

To find the lowest energy eigenfunctions and eigenvalues of the energy-density operator in 
(2) we substitute unknown functions Rn{z) in the form (12) into expression for H(z) in (11), 
and insert the resulting expression for n(z) into the equation for eigenfunctions of the operator 
H in (2). We arrive at the following system of non-linear algebraic equations for the lowest 
eigenvalues En of the operator H and coefficients aj, 6,: 
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ailloo + ^illoi = aiEo 

a2noo + 621101 = a2-Ei ^ ' 

b2Wi + 621111 + a2nio = 62-Bi 

where the matrix elements Hmn of the "potential energy" operator 7rn(z) arc taken on the set 
of eigenfunctions Yn{z) defined above. It proves to be that when IIoo/M^i < 10 the system of 
equations (13) can be solved in closed form. The solutions are: 

Eo = noo(l - hS), Ei^Wi+ noo(l + h + /3<5), (14) 
(ai;6i) = (l + 52)-i/2(i._j)^ (02:62) = (1 + 52) -^/'(<5;1). (15) 

Here, the following parameters have been introduced: 

Is = / Y^{z)dz « -0.4, h = L^ I Yt{z)dz « 1.8, 
Jo Jo 

5 = /3noo/(Wi + noo/4) « -0.17. (16) 
The matrix element IIoo obeys normalization condition: 

nPeff 



Hoo = Y [ ■^'n{z)dz = 
i-' Jo 



Effective surface tension in the bilayer,Pe//; balances entropic repulsion between hydrocarbon 
chains. Substituting (15) into (12) and subsequently (12) into (11) we find analytical expression 
for the lateral pressure profile: 

n(z) = ^(^p^ {(i^o(^) - 6YMf + {Y,{z) + 6Y,{z)f) . (17) 

Normalized lateral pressure profile (17) is plotted in Fig. 3. The system parameters are: chain 
full length L = 20 A, area per chain in all- trans conformation A = 20^^, chain bending rigidity 
Kf = ksTlp = 2.8 • 10~^ erg-cm, T = 3007sr, ks- Boltzmann constant, Ip ~ i/3 - chain 
persistence length, Pej/ w 150dyne/cm^. The calculated profile, shown in Fig. 3, qualitatively 
agrees with the results of the recent molecular dynamics simulations*,'^. The behaviour of the 
11(2;) curve can be understood comparing fluctuations of a lipid chain in bilayer and in empty 
space. In the bilayer the chain fluctuations are strongly suppressed by the neighboring chains 
with respect to a free lipid molecule. Suppression of chain fluctuation leads to the decrease 
of conformational entropy and, in turn, results in the enhancement of entropic; pressure. The 
more suppressed the fluctuations, the higher the entropic pressure. The most significant relative 
decrease of fluctuation amplitude takes place at chain free end (2; ~ i); at head group region 
{z ^ 0) the relative decrease is noticeable too as lipid heads are not fixed; at mid-chain region 
the decrease of fluctuations is less pronounced due to segment connectivity in the chain. Now we 
can validate the exploited approximation of small chain deviations in the bilayer. We calculate 
the thermodynamic average of chain fluctuation amplitude using expansion (3) and distribution 
function (8): 

< R\z) >= 2 < Rljz) >= UbTY, (18) 

n 

Since, cx n'^, the sum in (18) converges fast. Allowing for relations Ro{z) ~ 1/L,Eq ^ IIoo ^ 
Peff/L, we find rough estimates for the small parameter: 

V< RHz) >/L ~ {kBT/L'^P^ffYI^ ~ 0.1. (19) 

To summarize, we have demonstrated that lateral pressure profile in the hydrophobic core 
of lipid bilayer is expressed via the eigenfunctions of the density operator of the chain energy. 
The self-consistent solution for the lateral pressure distribution is found. 
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FIG. 1. Mean-field approximation for inter-chain interactions in lipid bilayer. 
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FIG. 2. Hydrocarbon chain as a flexible string of finite length. R is vectorin {x, y} characterizing 
deviation of the chain from z-axis; A = n < >= nd^ is average area per string. 
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FIG. 3. Calculated lateral pressure distribution in the hydrophobic core of the bilayer. The pressure 
arises from entropic repulsion between hydrocarbon chains in each monolayer, z is coordinate along 
chain axis spanning from z = at head group to z = L at free chain end. Pressure is normalized by 

Ho = Peff/2L. 
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